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ABSTRACT 

Constructing dynamical maps from the filtered output of numerical integrations, we 
p I ■ analyze the structure of the vq secular resonance for fictitious irregular satellites in 

r-pl | retrograde orbits. This commensurability is associated to the secular angle 9 — zu~zuq, 

where w is the longitude of pericenter of the satellite and Wq corresponds to the (fixed) 
'-^ ■ planetocentric orbit of the Sun. Our study is performed in the restricted three-body 

problem, where the satellites are considered as massless particles around a massive 
Q . planet and perturbed by the Sun. 

S-^ ' Depending on the initial conditions, the resonance presents a diversity of possi- 

, ble resonant modes, including librations of 9 around zero (as found for Sinope and 

Pasiphac) or 180 degrees, as well as asymmetric librations (e.g. Narvi). Symmetric 
modes are present in all giant planets, although each regime appears restricted to 
certain values of the satellite inclination. Asymmetric solutions, on the other hand, 
seem absent around Neptune due to its almost circular heliocentric orbit, 
f — | Simulating the effects of a smooth orbital migration on the satellite, we find that 

O^l i the resonance lock is preserved as long as the induced change in semimajor axis is 

much slower compared to the period of the resonant angle (adiabatic limit). However, 
the librational mode may vary during the process, switching between symmetric and 
asymmetric oscillations. 

Finally, we present a simple scaling transformation that allows to estimate the 
resonant structure around any giant planet from the results calculated around a single 
primary mass. 
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1 INTRODUCTION In the same year, Saha & Tremaine (1993) found that 

the Jovian moon Sinope also displays a similar resonance 
Both Jupiter and Saturn are rich in irregular satellites. , , . 1x , . . , 

° lock m the resonance, altough m this case the motion al- 

Jupiter has at least 54 moons m this group, while Saturn , . i ■ i i r 

, ternates between hbration and circulation on timescales or 

contains more than 30 members. A smaller population is , r lr .s .„r „ , n . , n , .. , n 

the order or 10 — 10 years. Extending the orbital evolu- 
also observed around Uranus and Neptune, although it is r r lri 7 , , , / ri „ rioN , , , 
, . tion lor 5 x 10 years, INesvorny et al. (2003) showed that 



Pasiphae's libration is also temporary, and will switch to a 
circulation of the resonant angle in approximately 3 x 10 7 
years. 



not clear whether this is intrinsic or due to observational 
bias. 

Irregular satellites suffer strong perturbations from the 
Sun which significantly affect their orbits evolution. Whip- 
ple & Shelus (1993) numerically integrated the orbit of the Cuk & Burns ( 2004 ) Performed a detailed search for 
Jovian moon Pasiphae for 10 5 years under the effects of other irregular moons in secular resonances, finding that 
all planets, and found that its longitude of pericenter zu is the Prograde Saturn moon Siarnaq appears to the trapped 
locked in a secular resonance where vj-vj Tup librates around in a Pericenter secular resonance with Saturn (i.e. u 6 ), also 
180°. Here ot Jup is the longitude of pericenter of the helio- displaying intermittent libration, where the resonant angle 
centric orbit of Jupiter. Since vj Jup = vj q + tt, where vj & oscillates around 180°. While the Uranus irregular satellite 
is the Jupiter-centric longitude of the pericenter of the Sun, Stephano is not in libration, the angle w - ro shows a 
this implies that zu — zuq oscillates around zero. By anal- quasi-resonant behavior with a very long-period circulation, 
ogy from asteroidal motion, we will call this the V5 secular Finally, Beauge & Nesvorny (2007) found that although 
resonance for satellite orbits. the retrograde Saturn moon Narvi is not currently in reso- 
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nant motion, long-term simulations show a future libration 
in the vq secular resonance. Contrary to previous cases, here 
the resonant angle displays libration around values different 
from zero or 180°, in what appears to be an asymmetric 
libration point. Since all these resonant configurations, in- 
dependently of the planetary mass, involve the relative be- 
havior of the longitude of pericenter of the satellite and that 
of the planetocentric orbit of the Sun, in all future references 
we will denote this commensurability as vq . 

At present it is not clear whether the observed resonant 
population in the outer planets is evidence of a past smooth 
orbital migration of the satellites, or whether it is simply 
due to chance. Beauge & Nesvorny (2007) found that sta- 
bility criteria alone yield satellite populations close to vq, 
although the proximity of some individual moons (particu- 
larly Pasiphae and Sinope) appear sufficiently detached from 
a random distribution to be statistically significant. 

One of the problems in understanding a possible rela- 
tionship between the z/q resonance and the past evolution of 
the irregular satellites, is the lack of a detailed analysis of the 
resonant structure. As will be shown in Section 2, analytical 
models, even of high-order, are not sufficiently precise and 
are unable to reproduce the main characteristics of the sec- 
ular commensurability. Semi-analytical models, where the 
secular perturbations are evaluated by numerical averaging 
of the exact Hamiltonian, also suffer the same limitations 
since they are equivalent to a first-order theory. For this 
reason, so far most dynamical studies have been restricted 
to numerical integrations of individual orbits. 

In this paper we present a numerical study of the struc- 
ture of the vq resonance for retrograde orbits. The resonant 
behavior is obtained applying two filters on the numerical 
output: the first is constructed to eliminate the short-period 
variations (frequencies comparable with the mean motions). 
The resulting data is then filtered once again to reduce the 
effects caused by the non-resonant secular angular variable. 
The results are displayed as dynamical maps showing the 
location, type and extension of the librational domains for 
several values of the integrals of motion. 

The details of the numerical method are outlined in 
Section 3, which includes an application to the region of the 
phase space in the vicinity of the Jovian moon Sinope. Sec- 
tion 4 shows the structure of the vq resonance in the Saturn 
system. In both cases the maps are constructed for a single 
value of the proper semimajor axis. Section 5 discussed the 
effects of an ad-hoc migration acting of the satellites orbits, 
and its effects on the resonant motion. Next, in Section 6 
we present a simple scaling law that allows us to relate the 
resonant structure for any planetary mass. Finally, the ap- 
plication the real satellites are briefly discussed in Section 
7, while concluding remarks close the paper in Section 8. 



2 ANALYTICAL MODELS FOR THE v Q 
RESONANCE 

Suppose a fictitious satellite in orbit around a giant planet 
trip perturbed by the Sun. Let a be the satellite's planeto- 
centric semimayor axis, e its eccentricity, i its inclination 
with respect to the Laplace plane, M the mean anomaly, u> 
the argument of pericenter and Q the longitude of the node. 
Planetocentric orbital elements of the Sun will be identified 



by the index 0. Since we will adopt a Hamiltonian approach, 
it is usefull to introduce the modified Delaunay canonical 
variables: 

L = yfjm ; M 

G = LV1~^ ; oj (I) 
H = Gcosi ; fi. 

The complete Hamiltonian function can then be written as 
F = Fo(L, A) + R, where Fo is the two-body contribu- 
tion (A is the canonical momenta associated to Mq) and 
R — R(L,G, H,M,ui,£l, Mi) denotes the disturbing func- 
tion stemming from the Solar gravitational effects. R also 
depends on the orbital elements of the Sun (e©, zq, zuq 
and £2q), which are considered constant, as well as on the 
planetary (i.e. central) mass m p . 

Yokoyama et al. (2003) presented an analytical model 
for the secular behavior of irregular satellites in the re- 
stricted three-body problem, and applied it to the Jovian 
satellite system. It is based on Kaula's (1962) expansion 
of the disturbing function, truncated at fourth order in 
a = a/a Q and eccentricities and inclinations. The averaging 
over the short-period terms, associated to the mean anoma- 
lies of both bodies, is done to first order. We will denote as 
F\ the resulting expression of the secular Hamiltonian. 

Although we will not write the expression ex- 
plicitly, functionally it can be written as Fi — 
fi(G* ,H* ,uj* , fi*;a*), where the new (starred) variables 
are the mean elements, which must not be mistaken with 
the original osculating variables, the quantity a* is related 
to the "proper" Delaunay momenta L* by the equation 
a* = L* 2 /fia Q . 

To study a given secular resonance, we must perform a 
canonical transformation to resonant variables. The vq com- 
mensurability is defined by the angle 9 = vj — vjq ; however, 
the relationship between the argument of pericenter uj and 
the longitude of pericenter zu is different when considering 
direct or retrograde orbits. Since most of the real irregular 
satellites in the vicinity of this resonance are retrograde, we 
will assume i > 90°, for which zu — £2 — uj. 

Our set of resonant variables will be (ii, h, 4>, Q), where 
(f> = fi* — uj* will be a slow angle related to 9, while Q is 
a fast angle. This can be chosen among any linear combi- 
nation of the secular angles with the only condition being 
that it be independent of 4>. For the present work we will 
choose Q — Q* , and check that the precessional frequency 
of the longitude of the node is higher than the frequency 
of <f>. Having specified the angles, the expressions for the 
corresponding canonical momenta can be easily found. The 
complete transformation from the mean Delaunay variables 
to the resonant counterparts are given by: 

h = -G* ; 4, = Q?-u* 

h = G*+H* ■ Q = Q.*. {> 

The resonant angle would then be 9 — cf> — tuq . 

The transformation of the secular Hamiltonian to res- 
onant variables is straightforward, yielding F\(Ii, I2, 9, Q). 
Since Q is a fast variable, we can perform a first-order av- 
eraging of the Hamiltonian over this angle. This averaging 
implies a new canonical transformation to new resonant vari- 
ables (7j* , I2, 9* , Q*) such that the transformed Hamiltonian 
Fj — F*(Ii , 9*; I2) and Q* is cyclic. The associated mo- 
menta I2 is then a new constant of motion, and the system 
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is reduced to a single-degree of freedom model for the vq 
resonance in the canonical pair (/*,$*). 

Equation (10) of Yokoyama et al. (2003) gives explicit 
expressions for F*, that can be succinctly written as: 



Fx = A + Ax cos (6*) + A 2 cos (26»* 



(3) 



where the coefficients Ak = Ak(I*; L* , 1%) are function of 
the resonant momenta I{ and parametrized by the constants 
L*,l2 and the orbital elements of the Sun. Different initial 
conditions will give different values of the proper elements 



2.1 The Secular Dynamics of Sinope 

As an example, Figure [Ha) shows the level curves of con- 
stant F* in the plane (9* , e*), adopting Jupiter as the central 
body (i.e. m p ). The quantity e* = e*(I*; L* , 7|) is the so- 
called mean-mean eccentricity, which can be obtained from 
the canonical momenta I* assuming that the values of L*, 1% 
are constant for all initial conditions. The adopted values for 
these parameters are close to the irregular satellite Sinope, 
and were calculated from a numerical simulation of the ex- 
act equations (three-body problem) over a timespan of 10 5 
years. 

For these constants, the phase plane of F* shows a res- 
onance region centered approximately at e* ~ 0.85 corre- 
sponding to a libration of 9* around 180°. The width of 
the libration zone is approximately Ae ~ 0.05. All initial 
conditions with lower values of the eccentricity correspond 
to circulations. This structure contrasts significantly with 
the Figure die) which shows the real dynamical behavior 
of Sinope obtained with the same numerical simulation de- 
scribed previously. The values of (9* , e* ) were calculated us- 
ing a low pass FIR digital filter (Carpino et al. 1987) to 
eliminate all periodic variations with period smaller than 
100 years. This includes variations in the orbital elements 
due to both the mean anomalies and f2. For comparison, in 
gray dots we also show the evolution of the osculating ele- 
ments. A comparison between both figures shows that the 
analytical model F* fails to reproduce the correct dynamics 
of the system, predicting that Sinope should be far from the 
resonance domain displaying a circulation of 9*. 

This imprecision is not restricted to the structure of 
the phase space in vicinity of this satellite, but is a general 
flaw of the averaging process used to obtain the analytical 
model. As shown by Cuk & Burns (2004), the averaging over 
short period terms must be performed to order higher than 
unity, to include the effect of second and third order sec- 
ular terms such as evection. These effects cause significant 
variations in the secular frequencies, thus affecting both the 
long-terms evolution of individual bodies, calculations of in- 
cluding proper elements, as well as the location of secular 
resonances in the domain of irregular satellites. 

With this in mind, Beauge et al. (2006) developed a 
high order analytical model of the secular Hamiltonian for 
irregular satellites, where the averaging process over the 
mean anomalies was extended to third order. The new 
Hamiltonian function, which we will refer to as F3 = 
Fs(G* , H* ,uj* , Q*; a*) was used to calculate the proper ele- 
ments of all irregular satellites of the outer planets and esti- 
mate their proximity to different secular resonances (Beauge 
& Nesvorny 2007). 
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Figure 1. (a). Level curves of constant Hamiltonian F? for the 
Jovian vq resonance, in the plane (e*,6*), using a first-order av- 
eraging over short-period terms. The adopted values for the inte- 
grals of motion are L* = 2.115 X 10 _4 ,i"| = 1.575 X 10 -5 , in units 
of Solar mass, AU and day. These values correspond to the Jo- 
vian satellite Sinope. (b). Same as before, but for the third-order 
Hamiltonian F£ . (c). Numerical simulation of Sinope (three-body 
problem) over 10 5 years. Gray dots show the evolution of the oscu- 
lating elements, while black dots correspond to an output filtered 
to eliminate all periodic variations with period smaller than 100 
years. 



From the analytical expression of F3, we can now con- 
struct a higher order model for the vq resonance in the same 
manner as before, and plot the level curves of F£ in the vicin- 
ity of Sinope. Results are shown in Figure [lib) and show a 
significant improvement with respect to Fx- The resonance 
domain is now much closer to the actual orbit of the Jovian 
satellite, and the main libration island is now (correctly) 
associated to an oscillation around 9* — 0. However, some 
important discrepancies still remain. The center of the li- 
bration region is still not correct, and the analytical model 
now predicts a smaller eccentricity than the numerical sim- 
ulation. More important, _F 3 * shows an secondary libration 
island around 9* = 180° which in fact does not exist. As 
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Figure 2. Frequencies of the secular angles ut* , Q* and <j> = 
f2* — lu* , obtained with both analytical Hamiltonians (continuous 
and dashed curves). Open circles show numerical results calcu- 
lated from a Fourier analysis of filtered data. All initial conditions 
correspond to the same values of L* and /| as in figure [T] 



shown in Figure 16 of Beauge & Nesvorny (2007) the size of 
this second libration mode increases for larger values of 7|, 
becoming an important characteristic of the phase space for 
librations at high eccentricities. 

2.2 Estimation of the Secular Frequencies 

The top and middles plots of Figure [2] shows the frequencies 
of the secular angles ui*, Q* (denoted by v w and un, re- 
spectively) as a function of the proper eccentricity e*, with 
initial conditions corresponding to the same values of the in- 
tegrals L* and 7| as before. Dashed curves show the results 
obtained using the analytical Hamiltonian F* restricted to 
those terms which depend only on the canonical momenta. 
Frequencies calculated with the third-order function F£ are 
shown in continuous lines. For comparison, we also estimated 
the same frequencies numerically (open circles), Fourier an- 
alyzing the filtered output of a series of simulations for 10 5 
years. 

For both ui* and Q* there is a very marked offset be- 
tween the real frequencies and the results obtained through 
F* , and the analytical estimates systematically yield lower 
values. A better agreement is noted using F 3 * , especially for 
the precessional frequency of the node, although a larger 
discrepancy is noted for the argument of pericenter. It is 
also important to recall that both analytical models were 



developed using Kaula's expansion of the disturbing func- 
tion, whose convergence is only guaranteed for eccentricities 
below ~ 0.6 (e.g. Wintner 1941). This seems to be the cause 
of the increasing imprecision of F£ for values of e* approach- 
ing this limit, but also implies that the results obtained with 
the simpler model F* are also undependable for e* > 0.6. 

The bottom graph of Figure [2] shows the values of 
vq. — Vu for the same initial conditions. Zero values (dot- 
ted line) correspond to the exact vq resonance, since we are 
considering constant longitude of pericenter for Jupiter's or- 
bit, this plot shows how the imprecision in the secular fre- 
quencies affect the resonance location. For F* the exact res- 
onance occurs for almost parabolic orbits, well beyond the 
convergence limit of the expansion of the disturbing func- 
tion, and thus its very existence appears unreliable. Once 
again better agreement is noted using Fi, although a sig- 
nificant error is still in evidence. Notwithstanding, since the 
continuous curve is very shallow for low eccentricities, even 
small differences in the model can cause important displace- 
ments in the location of the exact resonance. 



3 THE NUMERICAL MODEL 

Outside secular resonances, analytical models constitute ad- 
equate tools for the long-term dynamics of irregular satel- 
lites, and can be used to obtain fairly precise proper elements 
for real bodies (Beauge & Nesvorny 2007). They are also 
able to reproduce the structure and location of the Lidov- 
Kozai resonance (Lidov 1961, Kozai 1962) with good preci- 
sion. However, even high order analytical models seem insuf- 
ficient for the task of mapping the structure of other secular 
resonances, such as the vq commensurability. Not only may 
the location of the resonance domain may be imprecise, but 
the libration centers and secondary modes predicted by the 
models may be unreliable. 

In this section we will present a series of dynamical 
maps of the Jovian vq resonance, obtained numerically by 
filtering the outputs of a series of long-term simulations of 
fictitious particles. Each map consists of the dynamical evo- 
lution of a set of initial conditions with predefined values of 
the integrals (L* , 1% )• The main advantage of this procedure 
over analytical expansions is that we are no longer restricted 
to low to moderate eccentricities, and there are no approxi- 
mations in the modeling of the gravitational interactions. As 
before we will work in the realm of the restricted three-body 
problem and assume a fixed orbit for the planetocentric mo- 
tion of the Sun. 

For each map we begin with a set of initial conditions 
in osculating orbital elements for the particle plus the mean 
anomaly of the Sun (a, e, i, M, Mq , u), Q). Each orbit is in- 
tegrated for a time span of 10 5 years using a Burlisch-Stoer 
based N-body code. The output is subsequently converted 
to canonical resonant variables (L, I\ , I2, M, Mq, <j>, Q) and 
filtered to eliminate the variations associated to the mean 
anomalies M, Mq , as well those corresponding to the fast 
variable Q. Although the orbital period of both bodies are 
well known beforehand, the frequency of Q may depend with 
the initial condition and must be monitored closely. Thus, 
the decimation and length of the filter may change in each 
simulation. However, we have found that a cutoff period of 
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200 years and a filter of 1600 data points gave very precise 
results in most cases. 

The output of the filtering is a numerical approxima- 
tion to the transformed resonant variables (/j*, </>*; L* , 7J), 
analogous to those obtained with an analytical model, with 
one important difference. Although the digital filtering can 
be considered analogous to an averaging of the short-period 
angles, it does not reduce the number of degrees of freedom 
of the system. Consequently, any chaoticity in the original 
solution will be preserved in the filtered data. This implies 
that a regular (non-chaotic) solution will give a one dimen- 
sional curve in the (I*,(f>*) plane, while a stochastic trajec- 
tory will cover a two-dimensional region. 

To construct the map, each set of initial conditions in 
osculating variables must yield the same values of (L* , 7|) = 



L* = 



From perturbation theory we know that 
d 



I 



L+ dM X 
J 2 + 7^X : 



: fi(a,e,i,M,Mo,<j,n) 
f 2 (a,e,i,M,M e ,Lj,n), 



(4) 



where x(I J ,Ii,l2,M,Mo,4>,Q) is the generating function 
written in osculating variables evaluated at the initial condi- 
tions. Analytical models allow explicit construction of func- 
tion x, which can be inverted to give the necessary values 
of the osculating elements. Since numerical simulations do 
not give fi, the inversion of expressions Q must be done by 
iterations. 

The idea is to introduce small deviations in the osculat- 
ing variables a —> a + Aa and i — > i + Ai and leave all the 
other initial conditions intact. Introducing these expressions 
in and expanding /, in a first-order Taylor around the 
original values, we obtain: 

— Aa + -4- At 



L* + AL* 



L* + 



7 2 * + A7 2 * = 7 2 * + 



da 



di 



(5) 



Aa + ^Ai, 
aa oi 



where the partial derivatives of fi are evaluated at the initial 
orbital set. Equating L* + AL* = L* and 7 2 + A7| = 1%, 
equations ((5]l can be inverted to give the necessary values 
of Aa and Ai. The partial derivatives can be estimated nu- 
merically using finite differences, calculating how the inte- 
grals vary with very small deviations in each of the orbital 
elements around the original values. This procedure is not 
very precise, but can be used a successive approximations 
to give the desired solution. In general, a few iterations are 
sufficient to give values with relative precision of the order 
of AL*/L*, A/J/7J ~ 10" 4 . 
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Figure 3. Filtered output, in the (e* ,0*) plane, of a series of 
numerical integrations with initial conditions leading to a* = 0.16 
AU (L* = 2.115 X 10~ 4 ) and four different values of 7 2 (in units 
of 10 — 5 ). Each plot is roughly centered on the libration domain 
of the Jovian vq resonance. Chaotic motion is noticeable in the 
lower right-hand plot, especially for values of e* larger than the 
center of libration. 
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Figure 4. Values of i* and e* corresponding to the center of the 



libration width in each orbital element. 



3.1 Resonance Maps Around Sinope 

Figure [3] shows several dynamical maps constructed for 
L* — 2.115 x 10" 4 (corresponding to a* = 0.16 AU) and 
four different values of 1% . For simplification purposes, all 
numerical values of J| will be given in units of 10~ 5 . The 
lower left-hand frame roughly corresponds the initial con- 
ditions of Sinope. Instead of plotting I* as function of the 
resonant angle 8*, we converted the canonical momenta to 
equivalent mean-mean eccentricities e* and inclinations i* 
according to the relationships: 

i? = —G* = -LVl-e* 2 (6) 



J| = G*+H* =77(1 + cos i*). 

In all the plots we notice a main libration island centered 
around 6* = 0, altough the map for 7J = 1.0 also shows 
a small libration region centered at 6* = 180°. The values 
of i* and e* of the center of each main libration island are 
shown Figure [4] where the error bars indicate the libration 
width in each orbital element. The libration zone appears 
extensive in the eccentricity, although it decreases in size 
for larger values of 7J. However, the corresponding size in 
i* is very small, indicating that the v® commensurability 
appears very restricted in the inclination. 

Although most orbits in Figure [3] appear regular, the 
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Figure 5. Averaged MEGNO indicator (Y) for 200 initial condi- 
tions with 8* = in the four dynamical maps shown previously. 
Regular orbits are characterized by (Y) ^ 2 while larger values 
are indicative of chaotic motion. Unstable orbits were given a 
value of (Y) = 10. In each plot the libration region is shown in 
shade. 
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Figure 6. Structure of the uq resonance for initial conditions 
around the Saturn irregular satellite Narvi. Each plot shows fil- 
tered output, in the (e*,0*) plane, of a series of numerical inte- 
grations with initial conditions leading to a* = 0.13 AU and four 
different values of 7| (in units of 10~ 5 ). 



map for 7J = 3.0 shows some irregular trajectories which 
appear to indicate chaotic motion. To test this, we inte- 
grated 200 initial conditions for every value of 7| • Each had 
a different initial value of 7^ and the initial value of the reso- 
nant angle 8* was taken equal to zero. Total integration time 
was 10 6 years and, together with the equations of motion, 
we also calculated the MEGNO chaoticity indicator (Cin- 
cotta & Simo 2000). Figure [5] shows the values of averaged 
MEGNO number (Y) as a function of the mean-mean ec- 
centricity for each map. Regular orbits yield (V) ^ 2 while 
larger values are indicative of chaotic motion (see Cincotta 
and Simo 2000 for more details). Unstable orbits, leading to 
an ejection of the particle in the integrated timespan, was 
assigned a value of (Y) — 10. The extent of the libration 
island is shown as a shaded rectangle. 

For low values of 7| the chaotic region is restricted to 
high eccentricities and far from the secular resonance region. 
Thus, all commensurate trajectories are indeed regular, at 
least in the time interval analyzed. For 7^ = 1.5 the stochas- 
ticity of the resonance separatrix is clearly visible. However, 
for 7| = 3.0 the chaos is much more extensive and covers 
the Jovian Uq resonance zone. This helps explain why the 
dynamical map constructed for this value of 7| showed ir- 
regular trajectories in the (e*,8*) plane. A simple analysis 
of the output in osculating elements shows that the stochas- 
ticity is caused by the 6/1 mean-motion resonance with the 
planetocentric orbit of the Sun. As shown recently by Hinse 
et al. (2009), the retrograde satellite region in the Jovian sys- 
tem is affected by several mean-motion commensurabilities. 
Even though they appear to be of high order, the interac- 
tion with the secular resonance gives origin to a noticeable 
chaotic behavior. 



4 RESONANCE MAPS IN THE SATURN 
SATELLITE SYSTEM 

Since the Jovian satellite system is tainted with chaos stem- 
ming from mean-motion resonances, we switched planets are 
analyzed the structure of the vq commensurability in Sat- 
urn. Beauge & Nesvorny (2007) found that although the 
irregular satellite Narvi is not currently in a resonant con- 
figuration, the effects of the other outer planets will cause 
temporary trapping in the uq resonance in the future. Con- 
trary to both Sinope and Pasiphae in the Jupiter system, 
the resonant behavior of Narvi appears to correspond to an 
asymmetric libration of 9*. 

Taking the orbital elements of Narvi as a staring point, 
we constructed resonance maps for Saturn in an analogous 
manner as those constructed for Jupiter. In this case the 
proper semimajor axis was a* = 0.13 AU (corresponding to 
L* = 1.04 x 10~ 4 ), and we varied 7| from 1.0 to 3.0 (in units 
of 10 ). For each value of 7| we searched for initial con- 
ditions with different values of 7*. The results were filtered 
to eliminate the short-period terms associated to the mean 
anomalies and the variations due to the precession of the 
longitude of the node. Results are shown in Figures [5] and [7] 
The first figure corresponds to lower values of 7J and show 
significant similarities with the maps constructed for Jupiter 
(Figure [3]). Libration at low eccentricities is associated to a 
libration around 6* = 0, although for 1% = 1.6 we note that 
the symmetric solution turns unstable and bifurcates into 
two islands of asymmetric libration. 

The asymmetric structure appears more evident in the 
top-left hand plot of Figure [7] (7| = 1.8). There are now 
two separatrix, an exterior homoclinic curve which stems 
from 8* — 180° and surrounds the whole resonant domain, 
plus an interior separatrix centered around 8* = that en- 
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160 I 1 1 r- 

l* = 1.0 

155 -tf 1. 



H 1-2 



150 ■ 



140 - 



135 ■ 



I— I 1.4 



H 1.9 



I— I 



0.1 



0.2 



0.3 

e* 



0.5 



Figure 8. For L* = 1.04 X 10" 4 (i.e. a* = 0.13 AU) in the Saturn 
system, plot shows values of i* and e* corresponding to the center 
of the main libration island for different values of IX. Error bars 
indicate the libration width in each orbital element. 



compasses the asymmetric libration centers. The order is 
inverted in the next plot, constructed for 7| = 2.0, and now 
the critical curve associated to 9* — appears exterior in the 
resonance domain. As 7| increases the equilibrium values of 
9* move towards 9* = 180° until a new bifurcation is noted 
to 7| ~ 2.8. From this point onwards the libration returns 
to a symmetric solution, this time around 9* — 180°. 

It is interesting to note the change in the topology of the 
resonance as a function of 7| and, consequently, as function 
of the value of e* associated to the center of the libration 
domain. Previous studies of the vq commensurability only 
showed a libration around 9* = and, in the case of Narvi, 
a possible asymmetric libration. The complete structure is 
far more complex. Figure [S] shows the values of i* and e* 
corresponding to the center of the main libration island for 
different values of 7| . The same general trend is observed as 



in Figure [4] where once again the resonance region appears 
restricted to mean-mean eccentricities below e* ~ 0.55 for 
high values of the orbital inclinations. 

Finally, no significant chaos was found for these initial 
conditions, indicating that the uq resonance appears more 
stable for Saturn than for Jupiter. At least in part, this is 
due to the larger value of a© which pushes the dynamically 
significant mean-motion resonances much further from the 
planet. 



5 SATELLITE MIGRATION 

So far we have only analyzed the resonant structure in Sat- 
urn for a given value of L* , corresponding to a* = 0.13 AU. 
Since the construction of the dynamical maps is very time 
consuming, instead of repeating the calculations for other 
semimajor axes we have chosen to study the effects of a 
slow (adiabatic) migration on the satellite orbit. The plan- 
ets were not affected by this migration and remained in fixed 
heliocentric orbits. 

To simulate a smooth orbital decay we used the same 
procedure as described in Beauge et al. (2006), which con- 
sists in including a Stokes drag-like exterior force in the 
equations of motion. This additional non-conservative per- 
turbation can be modified to introduce any desired change 
in both semimajor axis and eccentricity. For all the integra- 
tions we adopted values such that the migration introduced 
no change in either the eccentricity or inclination, but only 
a slow change in a* with a characteristic time scale equal to 
r a — 10 7 years. This value is much larger than the librational 
period of 8* , thus guaranteeing the adiabatic approximation. 

Figure [5] shows a typical example of an inward migra- 
tion, where the initial condition corresponds to a* = 0.13 
AU and 7j> = 1.5 and a moderate amplitude of libration 
around 9* — 0. The integration was continued until the final 
semimajor axis equaled a* ~ 0.06 AU (closer to the planet 
than the real irregular moons) , and the output was filtered in 
the same manner as before. The top graph shows the varia- 
tion of the resonant angle as function of the semimajor axis, 
and the orbital evolution caused by the migration occurs 
from right to left. Although the particle remains trapped in 
the secular resonance throughout the integration, for smaller 
a* it switches from an oscillation around 9* = to an asym- 
metric libration, and finally, to a motion around 9* = 180. 
The change in amplitude after each bifurcation is due to the 
discontinuity in the action at these points. The bottom plot 
shows the change in i* and e*. during the migration. The 
thick black curve was obtained applying an additional fil- 
ter in the numerical data to eliminated the variations with 
period similar to that of the resonant angle. 

The same procedure was repeated for initial conditions 
for other values of 7| and extended to an outward migration 
of the satellite orbits. The resulting families in the (i*,e*) 
plane are plotted Figure [10] Each thick continuous curve 
shows the filtered evolution of the resonant orbit during the 
orbital migration, and can be considered as the approximate 
location of the families of zero-amplitude solutions for dif- 
ferent values of a*. Empty rectangles show the location of 
the solutions for a* = 0.06 AU while diamonds correspond 
to a* — 0.19 AU. For larger semimajor axes the orbital be- 
havior showed strong chaotic motion leading to an escape of 
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Figure 10. Thick continuous curves show the families of cen- 
ter of libration as obtained from inward and outward migration 
starting from a* = 0.13 AU. The initial value of 7| of each family 
is marked close to the initial position (empty circle). The inward 
migration was stopped at a* ~ 0.06 AU and the corresponding so- 
lution marked with an empty square. The outward migration was 
continued until a* ~ 0.19 AU. Here the final solution is marked 
with a diamond. Dashed curves show the interpolated families for 
the three values of a* discussed above. 



Figure 9. Simulation of satellite orbital decay. Initial condi- 
tion corresponds to a moderate amplitude 9* = libration with 
I£ = 1.5 at a* = 0.13 AU. Migration only affects the satellite 
semimajor axis with r a = 10 7 years. In the lower plot the thick 
black curve is the filtered output over the resonant angle. 



the particles in short timescales. The initial conditions (i.e. 
center of libration for a* = 0.13 AU) are shown in empty 
circles. Dashed curves correspond to interpolated values cal- 
culated from the numerical solutions for all three semimajor 
axis. 

The area inside the "wedge" -shaped region gives a fair 
idea of the region covered by the vq resonance for retrograde 
orbits in the Saturn system, at least for the interval of a* 
analyzed. For high inclinations the libration domain occurs 
for low to moderate eccentricities, and the exact location 
depends strongly on a*. However, as i* — > 90°, the resonance 
region shrinks until for i* < 130° it is practically restricted 
to a single curve for all values of the semimajor axis. 



6 SCALING LAW FOR OTHER SATELLITE 
SYSTEMS 

The structure of the uq resonance shows the same dynami- 
cal structure when computed for fictitious retrograde satel- 
lites in the Jovian and the Saturn system. Although con- 
structed for different values planets and values of L* , Figures 
[3] (Jupiter) and[B] (Saturn) appear basically the same, even if 
the values of J| an d the e* of the librational centers are not 
equal. As a consequence, we wish to analyze whether there 
exists a simple transformation (or scaling) that allows us 
to deduce the main characteristics of the secular resonance 
around any planet, starting from the information obtained 
for Saturn. 

Switching planets introduces two main effects: (i) a 
change in the semimajor axis a© of the planetocentric or- 




Figure 11. Families of the zero-amplitude resonant solutions, 
as function of a*, for retrograde orbits in the Saturn satellite 
system and six values of /J. Black curves correspond to libration 
around 9* = 0, red to oscillations around asymmetric libration 
centers and blue to libration around 9* = 180°. Dashed lines 
show interpolated curves for a* = 0.06, a* = 0.13 and a* = 0.19 
AU. Thick green curve shows the interpolated family of solutions 
obtained for fictitious Jupiter satellites with a* = 0.16 AU (see 
Figures [3] and |4]|. 



bits of the Sun (to a new value (l'q), and (ii) a change in 
the mass m p of the planet (to a new value m' p ). The outer 
planets also differ in eccentricity and inclination, although 
these parameters can temporarily be considered fixed. 

Allowing for a change in ciq is fairly straightforward. 
Since the secular Hamiltonian (averaged over short-period 
terms) only depends on the ratio a/a© (e.g. Murray & Der- 
mott, 1999), the semimajor axis of the planet only appears 
as a normalizing factor in the dynamical system. Thus, it 
should be expected that the dynamical structure calculated 
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at a* in the first planet should correspond to the one located 
at a* — T a a* in the second primary, where 

r a = 5k. (7) 

a© 

A more complex task it trying to model the changes 
caused by differences in the planetary mass m p . However, 
we can use an approximation of the restricted circular three 
body problem (RC3BP) known as Hill's Problem (Hill, 1878, 
1886, see also Szebehely 1967, section 10.4). Hill's equations, 
originally constructed as the basis of a Lunar Theory, are a 
simplified version of the RC3BP in which the more massive 
primary (i.e. Sun) is considered to be infinitely far away but 
its gravitational effects are still felt in the system. In a ro- 
tating coordinate system with position vectors r centered in 
the smaller primary (i.e. planet), it is possible to introduce 
a scaling r — > /i 1//3 r which yields equations of motion that 
are independent of the masses (at least to smallest order in 
jl). We denote fi as the mass ratio between primaries. Con- 
sequently, the complete dynamical behavior of the system 
can be deduced for any given value of p, and extended to 
other masses just modifying the scaling parameter. 

Hill's approximation assumes very small mass ratios, 
motion close to the smaller primary and circular orbits be- 
tween both massive bodies. Although these conditions are 
not rigorously satisfied in our case, its results may still be 
used as a first approximation. We may then introduce a sec- 
ond scaling factor P m defined as 




which constitutes a new transformation in the coordinates 
due to changes in the planetary mass. Thus, given the com- 
plete change from one primary m p to another of mass m' p 
with different heliocentric semimajor axis, the original dy- 
namical structure obtained at a* should be reproduced at a 
new value given by 

a* = (r a r m ) o*. (9) 

As a first test, Figure[TT]once again presents the families 
of zero-amplitude resonant solutions in the Saturn system 
(continuous black, red and blue curves) for smooth varia- 
tions of L* . However, now colors indicate the type of libra- 
tion: black is used for motion around 9* = 0, red for asym- 
metric librations and blue for oscillations around 6* = 180° . 
As already shown in Figure [5] migration can cause changes 
in the mode of libration. The dashed curves show the inter- 
polated location of the centers for a* — 0.06, a* = 0.13 and 
a* — 0.19 AU. Thus, for any other value of the proper semi- 
major axis it is possible, at least qualitatively, to estimate 
the locus of zero-amplitude solutions in the (e*,i*) plane 
and the corresponding type of libration. 

As an example, let us consider the resonant maps con- 
structed in the Jovian system for a* = 0.16 (Figure [3J). The 
locus of resonant centers shown in Figure [4] can be inter- 
polated and is shown in Figure [11] as a thick continuous 
green curve. The individual values are marked by open cir- 
cles with the corresponding values of J| placed adjacently. 
Applying the scaling law ((9} we find that in the Saturn sys- 
tem this curve should correspond to a* — 0.197 AU, placing 
it barely above the upper dashed line. Not only does the lo- 
cation of the green curve shows a very good correlation with 



Jupiter (a*=0.105 AU) Saturn (a*=0.130 AU) 




Figure 12. Resonant maps, around all giant planets, for initial 
conditions leading to same location of the libration center (black 
circle) in the (e*,i*) plane. Values of the proper semimajor axes 
are given on the top of each frame. While the first three maps 
are very similar, the low eccentricity of Neptune prevents the 
appearance of the asymmetric domain. 



the prediction of this simple scaling law, but there is also 
an agreement between the predicted type of resonant mode 
and the dynamical maps shown in Figure [3] It is important 
to mention that the eccentricity adopted for each planetary 
body were equal to their individual present values. 

In conclusion, the libration centers (both location and 
type) of fictitious satellites around Jupiter can be deduced 
from their counterparts in the Saturn system, simply apply- 
ing the scaling law ([9j) to the original proper semimajor axis 
a* and analyzing the curve corresponding to the new a'* in 
Figure [11] However, we still need to analyze the structure 
of the resonance lobes in more detail and extend the study 
to other planetary masses. Results are shown in Figure [T21 
where we present four dynamical maps, each for fictitious 
satellites around a different giant planet. Starting from Sat- 
urn, the values of L* and adopted for the other frames 
were deduced from the scaling law. The value of was cho- 
sen close to the onset of the asymmetric librations, and thus 
constitutes a region particularly sensitive to changes in the 
dynamical system. 

If the planetary eccentricities were taken the same, the 
resonant structure around all four outer planets would be 
identical. However, since the value of is different for each 
primary, the resonant structure must reflect this change. The 
location of the libration center in e* and i* appears practi- 
cally invariant, with e* ~ 0.379 in all four plots. Moreover, 
Jupiter, Saturn and Uranus have very similar orbital eccen- 
tricities, and therefore the type of libration (i.e. asymmet- 
ric) is maintained and the resonant domain has similar sizes. 
However, since es a t > eijra > ej up , the libration lobe sur- 
rounding the asymmetric solution is larger in Saturn than 
in Uranus, while Jupiter presents the smallest region of the 
three. 
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Although the differences in topology of the resonant 
domain around the first three planets are not very signifi- 
cant, the map constructed for Neptune shows no asymmetric 
lobe and the libration occurs around 9* = 0. This, however, 
may be understood in terms of the planetary eccentricity 
eNcp = 0.008, much smaller than the other giant planets. 
Let us recall the first-order analytical model F£ presented 
by Yokoyama et al. (2003) and summarized in equation @- 
Writting the coefficients Ai explicitly in terms of the per- 
turber's eccentricity e , and retaining only smallest orders, 
we can express this function as: 

F* = Co + Cie Q cos (0*) + C 2 e Q cos (26*), (10) 

where the new coefficients C% = Ci(I*; L* , 7| ,m p ) are now 
independent of e©. Although we have mentioned that this 
model is not a reliable representation of the dynamics of the 
uq resonance, the imprecision is only due to the values of the 
coefficients d. The functional form of the resonant Hamil- 
tonian stems from the D'Alembert rules of any expansion of 
the perturbing potential, which specifies that, for any integer 
k, the terms associated to a fc-order harmonic of 9* must be 
proportional to e . Thus, any other Hamiltonian function, 
whatever its origin, will have the same form as l |10j l, 

Libration centers are specified by the conditions 
OF* /dl* = dF*/d6* = 0. In particular, the second con- 
dition explicitly reads 

sin (0*) ^Ci + 4C 2 e cos (6*)j = 0. (11) 

This equation has two trivial solutions (9* = and 0* = 
180°) which correspond to symmetric librations. However, 
if 

4|C 2 |e >\Ci\, (12) 

the contribution of the second harmonic in the resonant an- 
gle is larger than the first, and asymmetric librations are 
possible (see Beauge 1994 for a similar analysis for exterior 
mean-motion resonances in the asteroidal problem). Since 
C\ and Ci vary differently as function of L* and 1% , we can 
explain why asymmetric librations are found for some values 
of the inclination and not for others. However, condition (|12[) 
depends linearly on the planetary eccentricity. Thus, if e 
is sufficiently small, asymmetric solutions are not possible 
and all librations are of the symmetric type. This appears 
to be the case of the Neptune satellite system, at least for 
the values of L* , 1% analyzed. 

As a final test, Figure [13] shows a comparison of the 
resonant structure of two systems: the plot on the right re- 
produces the dynamics around the asymmetric libration lobe 
in Saturn (top-right frame of Figure [12)) . while the left-hand 
graph was obtained for a fictitious Neptune in which the 
planet's eccentricity was elevated to the same value as Sat- 
urn. As expected, the resonant structure is now equivalent, 
confirming that the absence of asymmetric librations in Fig- 
ure [12] is in fact due to the planet's quasi-circular orbit. 



7 CONCLUSIONS 

We have presented a numerical study of the structure of the 
vq secular resonance for retrograde satellites around the gi- 
ant planets. This commensurability is associated to a libra- 




Figure 13. Same as Figure [T2l comparing the resonant maps for 
Neptune (left) and Saturn (right). However, now the eccentricity 
of both planets were considered the same, and equal to the oscu- 
lating value of Saturn. The resonant structure in both frames is 
identical. 



tion of the secular angle 9 = vo — tcr©, where w — Q — u> 
is the longitude of pericenter of the irregular moon and ro 
the corresponding angle for the planetocentric orbit of the 
Sun. We have shown that the resonant domain presents a 
diversity of libration modes, including both symmetric and 
asymmetric librations. Each mode appears restricted to cer- 
tain intervals of the mean-mean inclination i* of the satellite. 
Librations around 9 — occur for large values of i* (typ- 
ically i* > 150°), while oscillations around 9 — 180° seem 
to dominate for i* < 140°. Asymmetric solutions appear at 
intermediate inclinations, although planets with very small 
eccentricities (such as Neptune) may not contain asymmet- 
ric solutions at all. Calculating the mean-mean inclination 
for real irregular satellites shows a good agreement with this 
predictions. Pasiphae (i* — 148°), Sinope (i* — 157°) both 
display librations around 9 = 0, while the resonant angle of 
Narvi (i* = 141°) is in an asymmetric mode. 

Simulating a smooth orbital migration of satellite orbits 
have shown that the resonant lock is maintained as long as 
the adiabatic condition is met; in other words, the charac- 
teristic timescale of the orbital variation is much longer than 
the libration period. We have used these results to map the 
location and type of libration centers as a function of the 
proper semimajor axis, thus giving a more global view of 
the region of the phase space dominated by the vq com- 
mensurability. 

Finally, we presented a simple scaling law that can be 
used to relate the resonant structure in any the satellite sys- 
tem of planetary mass. Since Jupiter, Saturn and Uranus 
have similar orbital eccentricities, it is expected that the 
structure of the secular resonance should also be very simi- 
lar, except for a scaling in the semimajor axis of the irregu- 
lar moons. Neptune, however, has an almost circular orbit, 
which inhibits the formation of asymmetric libration points. 
For satellites around this planet, it is expected that sym- 
metric solutions should dominate. 

From the scaling law, it can be shown that a change in 
the semimajor axis of the satellite is equivalent to a similar 
one on the planet. In other words, the same results would 
obtained if the heliocentric semimajor axis of the planet suf- 
fered a slow migration while the planetocentric orbit of the 
moon was kept fixed. Although we have only done a few sim- 
ulations with a fast (non-adiabatic) migration, it appears 
that the resonant configuration is not so easily preserved, 
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and any originally librating body was usually expelled from 
the resonance. However, additional simulations are required 
to confirm is this behavior is general. 

In this paper we have considered only the restricted 
three-body problem and, consequently, neglected the grav- 
itational effects of the other giant planets. This is indeed 
an approximation, since these additional perturbations may 
cause significant orbital variations (e.g. Carruba et al. 2004, 
Cuk & Burns 2004). In particular, the temporary nature of 
the vq resonance lock of real satellites is due to the gravita- 
tional effects of the other planetary bodies. In consequence, 
the results presented here must be viewed only as a start- 
ing ground for the complete secular behavior for retrograde 
irregular moons. Surely the full picture will be far more com- 
plex. 
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